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The equation of state (EOS) for partially ionized carbon, oxygen, and carbon-oxygen mixtures at temperatures 
3 x 10 5 K < T < 3 x 10 6 K is calculated over a wide range of densities, using the method of free energy mini- 
mization in the framework of the chemical picture of plasmas. The free energy model is an improved extension 
of our model previously developed for pure carbon (Phys. Rev. E 72, 046402). The internal partition functions 
of bound species are calculated by a self-consistent treatment of each ionization stage in the plasma environment 
taking into account pressure ionization. The long-range Coulomb interactions between ions and screening of 
the ions by free electrons are included using our previously published analytical model, recently improved, in 
particular for the case of mixtures. We also propose a simple but accurate method of calculation of the EOS of 
partially ionized binary mixtures based on detailed ionization balance calculations for pure substances. 

PACS numbers: 52.25.Kn, 05.70.Ce, 52.27.Gr 



I. INTRODUCTION 

An understanding of the physical properties of matter at 
high densities and temperatures is important as a problem 
of fundamental physics as well as for various applications. 
A particularly challenging problem is the calculation of the 
equation of state (EOS) for stellar partial ionization zones, 
where the electrons and the different ionic species cannot be 
regarded as simple mixtures of ideal gases: Coulomb interac- 
tions, bound-state level shifts, pressure ionization, and elec- 
tron degeneracy must be taken into account. In a previous 
publication [1], we calculated the EOS for carbon at temper- 
atures 10 5 K < T < 10 K over a wide range of densities 
p, based on the free energy minimization method [2, 3] which 
enables us to include the complex physics in the model and 
ensures thermodynamic consistency. For the case of carbon, 
in particular, the EOS developed by Fontaine et al. [4] (FGV) 
more than three decades ago and widely used in astrophysics 
up to now is based on the free energy minimization method at 
relatively low densities p < 0.1 g cm -3 . 

The free energy model inevitably becomes complicated 
when density increases above p > 0.1 g cm -3 because of the 
growing importance of nonideal contributions and the onset 
of pressure ionization. This latter phenomenon is difficult to 
treat in the framework of the "chemical picture" of plasmas, 
which assumes that the different ion species retain their iden- 
tity (see, e.g., Refs. |5|-|7|], for discussions). On the other hand, 
EOS calculations within the more rigorous "physical picture", 
quite successful at relatively low p (e.g., become pro- 
hibitively complicated at high densities. First principle ap- 
proaches based on path integral Monte Carlo(PIMC) HI O 
or molecular dynamics (MD) calculations ifTTl [l2tl are com- 
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putationally highly expensive, especially at high temperatures 
when excited ionic cores should be considered. To the best of 
our knowledge, no PIMC or MD data for carbon or oxygen 
in the temperature-density range of interest in this paper has 
been published so far. 

Besides these first principle based computations and a few 
works that consider the influence of anisotropic distributions 
of neighboring ions (e.g. Refs. Il3l ll4|P . hot dense plasma 
studies of medium- to high-Z elements are generally built 
around the modeling of an ion in its plasma environment plac- 
ing it at the center of a spherically symmetric system. The 
complexity of the physics that is still required leads to a wealth 
of models. In ion sphere (IS) models lfl5ll neighboring ions 
act as a mere neutralizing background beyond some radius, 
while some correlation with the central ion can be introduced 
through the computation of pair distribution functions flal . 
for instance in the hypernetted chain approximation (HNC) 
Inll . The free electron background ma y b e obtained from 
model prescriptions like Thomas-Fermi 11181 Il9ll or may in- 
volve quantum computations, usually from the Kohn-Sham 
equations in the density functional theory (DFT) context 
II20[ ulh . Several levels of refinement are possible to model 
the atomic structure of the central ion. Most of them merge 
the various excitation and/or ionization states into a fictitious 
average atom (AA) l22ll . This saves a lot of computation 
time because the overall self-consistent scheme for the ion 
and its environment has to be solved only once for a given 
thermodynamics condition. The price to pay, especially when 
interested in opacities, is some additional statistical proce- 
dure to infer individual ionization stages populations from the 
average-atom solution (e.g. B23I1 ). This undoubtedly works 
well for high-Z elements due to their innumerable quantum 
states which translate into unresolved transition arrays (UTA), 
but for lighter elements such as carbon or oxygen the rele- 
vance of the AA scheme becomes questionable. 

In this paper we employ a generalized version of the EOS 
model [1] which relies on the free energy minimization in 
the framework of the chemical picture and is applied to high 
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densities across the pressure ionization region. We com- 
bine separate models for different ionization stages in the 
plasma environment, taking into account the detailed structure 
of bound states (configurations, LS terms), and use Boltz- 
mann statistics to sum up the internal partition functions of 
these ions. The atomic structure of each ion embedded in 
the dense plasma is calculated using a scheme [24] which 
self-consistently takes into account the modification of bound 
states due to the environment. The free electron density 
around each ion is treated quantum mechanically, thus re- 
solving the resonances. Though neighboring ions act at this 
level as a neutralizing background, the long-range interac- 
tions in the system of charged particles (ions and electrons) is 
included in the thermodynamics using the theory previously 
developed for fully ionized plasmas (see Ref. 02511 and refer- 
ences therein). This model allows us to obtain not only the 
thermodynamic functions, but also directly number fractions 
for every ionization stage, unlike A A models. 

As different ions are treated on an equal footing, whether 
from the same element or not, our approach allows to treat 
mixtures in a similar way as pure plasmas. We apply the 
model to carbon and oxygen plasmas at temperatures between 
3 x 10 5 K and 3 x 10 6 K and mass densities in the range 
10 -3 g cm -3 < p < 10 4 g cm~ 3 . We also consider carbon- 
oxygen mixtures at the same plasma parameter values and 
propose a simple but accurate method of calculation of ther- 
modynamic functions of the mixtures based on the solution of 
the ionization equilibrium problem for the pure substances. At 
T < 3 x 10 5 K the model remains valid in principle, but the 
calculation becomes numerically much more difficult in the 
current implementation, and it has not been realized in this 
work. For temperatures higher than ~ 3 x 10 6 K our model 
recovers a simpler one where the free electron density is as- 
sumed to be uniform. 

In Sec. [n] we briefly describe the total free energy model. 
In Sec. [Ill] we present the model for internal free energy of 
ions, which takes into account bound state configurations in 
LS coupling and their interactions with the continuum of free 
electrons. The technique for the calculation of thermody- 
namic functions at equilibrium and its update with respect to 
Ref. llJ] is briefly described in Sec. [IV] In Sec. [V] we dis- 
cuss the results of the EOS calculations for carbon and oxy- 
gen plasmas and for their mixtures. Section [VI] is devoted to 
the conclusion. 



II. FREE ENERGY MODEL 



written as 



(1) 



where Xj = (2nft •? /mjk^T) \ ' is the thermal de Broglie 
wavelength of the ions of the jth chemical element in the 
plasma, mj is the mass of these ions, and fee is Boltzmann 
constant. For the electrons at arbitrary degeneracy, F e can be 
expressed through Fermi-Dirac integrals (we calculate F e and 
its derivatives using the code described in B25I1 ). The nonideal 
term can be written as 
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where the first term, Fc ou h includes contributions due to the 
long-range part of the Coulomb interactions between different 
(classical) ions, between free electrons (including exchange 
and correlations), and between ions and free electrons. The 
second term, Fi nt , linked to internal partition functions, in- 
volves sums over localized bound states around the nuclei and 
includes interactions between bound and free electrons. No 
strict definition of either free and bound electrons or ions ex- 
ists in a dense plasma; therefore, in general, the terms in Eq. 
d2j are interdependent. In our approach, we handle this diffi- 
culty as follows: we first calculate the properties of the ions, 
treated individually, embedded in the plasma, by developing 
self-consistent models for these ions, we then couple these 
models with a model which describes the long range interac- 
tions, and finally we minimize the resulting total free energy 

We calculate the Coulomb term and its thermodynamic 
derivatives using previously published fitting formulae (see 
Ref. f25ll for references). Compared to our previous paper 
13], there are two main improvements in the calculation of 
this term. First, we employ a correction to the linear mix- 
ing rule for the ions of different types, derived in l26ll . Sec- 
ond, we have implemented fully analytical calculations of all 
derivatives of Fc ou i needed to obtain thermodynamic func- 
tions (previously some derivations were made numerically). 
The latter improvement increases the accuracy of the calcu- 
lations of these derivatives and allows us to proceed straight- 
forwardly, with no need to include the additional refinement 
described Sec. IIIC of Ref. Ill], based on an extraction of 
the long-distance contribution -Fboui from calculated values 
ofF tot . 



Consider a plasma consisting of N e free electrons and Ni — 
Y^jv Njv heavy ions in a volume V, where Nj v is the number 
of ions of the jth chemical element having v bound electrons, 
and v can range from to Zj, where Zj is the jth element 
charge number. The free energy model is basically the same 
as in Ref. Ufl. The total Helmholtz free energy is F to t = F e + 
Fi + F cx , where Fi_ e denotes the ideal free energy of ions and 
free electrons, respectively, and F cx is the excess (nonideal) 
part, which arises from interactions. The term Fi is the kinetic 
free energy of an ideal Boltzmann gas mixture, which can be 



III. BOUND-STATE CONTRIBUTION TO THE FREE 
ENERGY 

In order to evaluate F[ nt , we calculate the ionic structure 
in the plasma using the scheme described in 12411 . It is based 
on the ion-sphere approximation, which replaces the actual 
plasma environment for every ion by the statistically averaged 
plasma effects on the electron wave functions within a spher- 
ical volume centered at the ionic nucleus. At present we do 
not include neutral atoms (v = Zj for the jth element), which 
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is justified at the temperatures and densities where the ioniza- 
tion degree of the plasma is high. For each ion containing v 
bound electrons, a radius of the ion sphere Rj V is determined 
self-consistently from the requirement that the sphere is over- 
all electrically neutral. The Hamiltonian describing the ion 
(j, v) immersed in the plasma is written as 

V 

H 3 v = ^h lv { ri ) + W ]v , (3) 

i=l 

where 

Vj v is the potential due to the plasma on the ion (j, v), that 

must be determined self-consistently, and V£" is a scaled 
Thomas-Fermi potential of the nucleus and v bound electrons 
ll27ll independent of the density and temperature. Note that 
Vl^ disappears in Eq. ([3). The effective Hamiltonian hj U gen- 
erates a one-electron wave function basis with a finite number 
of bound states. The coordinate parts ip 3 7 ^ lm of the bound-state 
wave functions are obtained from the Schrodinger equation 

h 3^nlm = e ^«^L, (6) 

where n,l,m are, respectively, the principal, orbital, and mag- 
netic quantum numbers for a given orbital. Then Hj V is diag- 
onalized in the subspace of Slater determinants generated by 
the set of all ^djL. The interaction term Wj V is responsi- 
ble for the LS splitting of configurations. It may also lead to 
some configuration interactions. As Wj V does not depend on 
the plasma properties, its matrix elements are influenced only 
through the modifications of the wave functions in Eq. (O due 
to plasma effects through Vj". Our experience shows that the 
i/-electron energies of the bound states of the ions are well 
approximated as E jva = E] va + E(„i)ea( e jW - e°- vnl ), 
where E® va and e^ vnl are calculated for the isolated ion, and 

a = (nl)i(nV)2 ■ ■ ■ {nl) u 2S+1 L defines a particular LS term 
of a configuration. The boundary condition at the ion sphere 
radius Rj V to solve Eq. (|6]l does not noticeably affect Ej ua 
except near the densities where the corresponding term a be- 
comes pressure-ionized. At these densities we use the depen- 
dence of the one-electron eigenenergy on the external bound- 
ary condition in order to estimate the degree of electron der- 
ealization due to pressure effects and to determine the corre- 
sponding occupation probabilities, in the same manner as in 
Ref. fll. 

The free electron density n^(r) and the potential Vj v (r) 
are determined self-consistently in the potential generated by 
the nucleus charge and a Boltzmann average of the ^-electron 
wave-functions, solving for Kohn-Sham states in the local 
density approximation of the DFT, and assuming Fermi dis- 
tribution of the free electrons at a given chemical potential 
jj, e . The contributions from resonances are taken into account, 



1 1 1 1 1 1 1 1 1 

" V 


1 | 1 , 1 1 , 1 1 , 1 | 1 1 , , | 1 1 1 , | , 





— ^^y' : 


" 1 




2 




3 


Carbon'05" 


- 4 


present ~ 


- 5 — ' - 
, 1 , , , , 1 , 


T = 2.2xl0 6 K 
, , 1 , , , , 1 , , , , 1 , , , , 1 , , , , 1 , 





21 22 23 24 25 26 27 



l°gio n e (cm" 3 ) 

FIG. 1: (Color online) Neutrality volumes v v normalized with the 
fiducial number density of free electrons n° e (see text), as functions of 
n° for T — 2.2 x 10 6 K for carbon ions with v bound electrons (v = 
0, 1, . . . , 5). The results from Ref. 1 1] (dashed lines) are compared 
with the present improved results (solid lines). 



as explained in Refs. [Qldil]. Together with v?? and Vj v ' , the 
ion sphere radii Rj V and the corresponding neutrality volumes 
Vj U = 4ttR^ v /3 are obtained from the neutrality condition for 
each ion. 

For a uniform free electron background, we would have 
v jv = v° jv = (Zj - v)/n°, where n°(/i e ,T) is the "fiducial" 
electron density, which would correspond to the true electron 
density of the ideal Fermi gas of electrons at the given \i e and 
T values. When taking into account the interactions of the 
free electrons with the ions, Vj V deviates from Vj u , as illus- 
trated in Figs.Q]and[2]for carbon and oxygen plasmas, respec- 
tively. The drops of the curves at certain values of n®([i e , T), 
which are especially sharp at lower temperatures (see Fig. [2]), 
are the consequence of pressure ionization of separate bound 
levels: when a nl level of ion crosses the continuum limit 
and appears as a resonance in the neighboring ionization state 
v, the corresponding ion sphere shrinks to compensate this in- 
crease in the free electron density of states (note that Z— v free 
electrons are still present in the volume Vj u ). For a given level 
nl, pressure ionization occurs at electronic densities which 
roughly depend on the ionic charge as (Zj—v) A . If the number 
of bound states that remain in the ion sphere is not sufficient 
to accept v electrons, that ion disappears: v = 5, 6, 7 disap- 
pear with the 2p sub-shell, and ^ = 3,4 with 2s, see Figs. Q] 
and [2] We emphasize that, in contrast to average-atom mod- 
els, separate ions have their own fate. Compared to Ref. [1], 
we have extended the calculation of the free electron states to 
higher electron energies and to higher /j, e . The consequences 
of these improvements for the neutrality volumes are seen in 
Fig. Q] where the dashed lines show our previous results and 
the solid lines correspond to our new calculations. The exten- 
sion of /i e to higher values allows us to capture the last dips 
of the neutrality volumes of naked nuclei {v — 0) as functions 
of rig, which appear at nP e > 10 26 cm -3 for both carbon and 
oxygen. These dips are due to the coupling of the hydrogen- 
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FIG. 2: (Color online) Neutrality volumes v u normalized with rv e 
(analogous to Fig.Q}, as functions of n° for T = 2.2 x 10 6 K (solid 
lines) and T = 5.52 x 10 K (dotted lines) for oxygen ions with v 
bound electrons (i> = 0, 1, . . . , 7). 



like {v = 1) bound states to the continuum, and are related to 
the pressure ionization of the K shell (where the curves v = 1 
in Figs.Q]and[2]terminate). 

Pressure ionization of each atomic shell is accompanied by 
a spreading of the electronic level into a band, eventually pass- 
ing into the continuum. We account for this effect by using 
occupation probabilities Wj va — Y[(nl)ea w 3 l/nl ' which are 
equal to the fraction of the band corresponding to the elec- 
tron configuration a that remains below the continuum level 
at given (T, p e ) values (see Ref. |l|] for details). 

In order to get a smooth pressure ionization and a smooth 
extension from calculated values of the neutrality volumes to 
higher densities, we replace vj u by v* v — Vj V Wj V + Vj v (l — 
Wj V ), where Wj U is the occupation probability Wj ua calcu- 
lated for the ground-state configuration a of the ion with v 
electrons. 

The separation of H u into parts (|4]i and allows us to 
capture the plasma effects in one-electron energies and wave 
functions through Eqs. and (O, while the ^-electron struc- 
ture - configuration energy mean shifts and LS splitting - re- 
sults from the contribution W v . In Ref. |lf] we calculated this 
contribution to binding energies using a modified version of 
the "Superstructure" code [27], as in Ref. l24ll . But our expe- 
rience has revealed that the use of the atomic code is unnec- 
essary for our purposes. This contribution is calculated once 
for the isolated ion and then combined with the one-electron 
energies €j Vn i (see above) at any density. Indeed the depen- 
dence of this contribution on the plasma conditions is weak (of 
the second order in perturbation) compared to the dependence 
of the energies €j m i- In the present work we have updated 
the corrections due to W v by using the detailed atomic energy 
database of the Opacity Project (OP) (Ref. 112811 and references 
therein). 

Having calculated energies Ej va and occupation probabil- 
ities Wj V a at given values (T, p e ), we obtain the contribution 
of the internal degrees of freedom of all the ions to the free 



energy as 

Zj 

-Pint = -Nik B T*^2 X! x i v h 1 - 2 ^, ( 7 ) 

j u=0 

where Z ju = (2S+1)(2£ + 1) exp(-E jua /k B T) 

is the internal partition function of the ion (j, v) in the plasma, 
and Xj V = Nju/Ni its number fraction . v Xj V = 1). Since 
we do not consider neutral atoms, the actual upper limit of v 
summation is Zj — 1. 

IV. THERMODYNAMIC EQUILIBRIUM 

We use a generalization of the free energy minimization 
method of Ref. [1] to the case of multicomponent mixtures. 
The minimum of F to t is sought at constant V, T, and Nj. 
The free parameters are the chemical potential p, e and the ion 
fractions Xj V , subject to constraints 

Zj - 1 Zj - 1 

Xj v > 0, ^2 x ]v =*}> N i / ! XjvVjv = V, (8) 

v— j u—1 

where Yj are fixed chemical element abundances (Yj ^ 0, 
J2j Yj = 1). The last condition in Eq. (JHJ reflects that V 
is constant. We ensure this condition by using the Lagrange 
multiplier method. The other conditions in Eq. (JHJ are im- 
posed explicitly in our minimization procedure by allowing 
the set of {xj U } to contain only the values for which these 
constraints are satisfied. Note that the charge neutrality condi- 
tion is fulfilled automatically, because the volume Vj V around 
every nucleus is neutral by construction. Some further details 
of the minimization algorithm can be found in Ref. [ 1 ] . For a 
mixture of chemical elements, we have in total l+X^j {%j ~ 1) 
independent parameters of minimization. 

The minimization procedure provides the values of F to t 
with a typical accuracy of four digits. Further improvement 
of the accuracy becomes problematic. The main sources of 
the numerical noise are the finite precision of calculation of 
the internal partition functions and neutrality volumes and the 
finite precision of the minimization procedure. The achieved 
accuracy is not sufficient for an accurate evaluation of thermo- 
dynamic functions, especially those that involve second and 
mixed derivatives of F to t (specific heat Cy, logarithmic pres- 
sure derivatives with respect to density, \ p , and temperature, 
Xt, adiabatic temperature gradient, and so on). To overcome 
this difficulty, we use an additional filtering of the numerical 
noise after the minimization. The filtering procedure is mod- 
ified relative to that in Ref. QJJ] . Specifically, for the results 
shown in the next section, we have performed calculations for 
a grid of points (p,T) separated by Alog 10 p = 0.05 and 
A log 10 T = 0.1, and at each point we have determined an im- 
proved value of Ftot using a bicubic 10-parameter polynomial 
fit to Ft t (p, T) on this grid. The fitting was done by x 2 min- 
imization with weights decreasing with increasing distance 
between a given point (p,T) and the grid point (p',T') in 
the log 10 p-log 10 T plane as [1 + (log 10 (p'/p)/A log 10 p) 2 + 



5 



(log 10 (T'/T)/A log 10 T) 2 ]" 1 . This improved method of pri- 
mary filtering of numerical noise allows us to avoid the addi- 
tional secondary smoothing of thermodynamic functions that 
was employed in Ref . |Q]] . 



V. RESULTS 
A. Ionization equilibrium 

The minimization of the total free energy immediately gives 
a solution to ionization equilibrium. Examples of partial frac- 
tions of ions as functions of T at constant p, or of p at con- 
stant T, are shown in Figs. [3]-|6] as well as in panels (b) of 
Figs. [7] and [8] in Sect. IV Bl Figure [3] shows current results 
(solid lines) in a pure carbon plasma for the partial fractions 
of ions with different numbers of electrons (marked by the 
numbers near the respective curves) and for the mean effective 
charge (marked by the symbol (Z) and scaled to the right ver- 
tical axis) at (a) p — 0.2 g cm -3 and (b) p — 10 g cm -3 . This 
effective charge is calculated by taking into account the par- 
tial derealization of the outer-shell electrons in the pressure 
ionization region as (Z) = ^2j V (Zj — v + 1 — Wj V ) Xj V . The 
ions recombine with decreasing temperature at both densities, 
but at the high density p = 10 g cm~ 3 there is no room for the 
existence in carbon of bound states other than Is (as can be 
inferred from Fig. Q] though for a fixed temperature: p = 10 
gem -3 corresponds to n e — 2 - 3 x 10 24 cm -3 .) Hence 
recombination can not proceed past the He-like (v = 2) stage. 

The dashed lines in the figure correspond to a simplified 
approach, where the neutrality volumes Vj U and internal parti- 
tion functions Zj„ are calculated assuming the uniform elec- 
tron density n e = n°. In this simple approximation we can 
obtain results in a wider range of T than in the accurate calcu- 
lation. This approximation works well at relatively low den- 
sities, as proved by the good agreement of the dashed and 
solid curves in Fig. [3ja). With increasing density, the cou- 
pling of the bound and continuum states becomes significant. 
More accurate calculation becomes necessary, as can be seen 
in Fig. Ob), especially at low temperatures where the differ- 
ence in the ion sphere between a uniform free electron density 
and the density constructed from quantum continuum states 
will be more pronounced. 

We also report on Fig. [3] the mean ion charge as obtained 
from FLYCHK (dash-dotted lines), a suite of codes based on 
simplified atomic models aimed at providing fast collisional- 
radiative models and spectra [29]. FLYCHK gives gener- 
ally a lower mean effective charge, especially at temperature- 
density conditions where several ionization stages are present. 
This is probably due to its treatment of density effects through 
a simplified approach. The experimental data points of Gre- 
gori et al. J30ll for (Z) are plotted in Fig.[3ja) as circles. They 
agree with our results at this density, especially if one consid- 
ers the experimental uncertainties of ±20 eV for the tempera- 
ture and ±0.25 for (Z) Bill . Experiments at higher densities 
would certainly be very valuable. 

The effects of pressure ionization are most appreciable in 
Fig.H] where ion population fractions (a) and the mean effec- 
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FIG. 3: (Color online) Number fractions of different carbon ions (left 
vertical axis, marked by the numbers of bound electrons v near the 
curves) and mean effective charge (Z) (right vertical axis) in a pure 
plasma versus temperature at (a) p — 0.2 g cm -3 and (b) p = 10 
gem -3 . Solid lines: accurate results; dashed lines: approximation 
where partition functions Zj V and neutrality volumes fj„ are cal- 
culated assuming a uniform free electron density. Dash-dotted line: 
results for (Z) from FLYCHK 1 29] . In panel (a) the circles are the 
experimental results from Gregori et al. [30] for p — 0.2 g cm -3 . 



tive charge (b) for a pure oxygen plasma (solid lines) are plot- 
ted versus density at T = 6 x 10 5 K. The results of FLYCHK 
for (Z) are also reported on Fig. @fb) (dot-dashed line). Re- 
combination with increasing density is modulated at the low- 
est densities (p < 0.5 g cm~ 3 ) by the successive pressure ion- 
izations of the Li-like (v = 3) ion shells (see Fig. [2 p = 0.5 
g cm -3 corresponds to n e ~ 10 23 cm -3 for (Z) — 5.5). The 
temperature is sufficient to populate excited states of that ion 
but not of He-like ions, so that only the Li-like ion partition 
function is affected (the Is 2 2s-ls 2 3s energy difference is 
9 x 10 5 K, while the ls 2 -ls 2s one is 6.5 x 10 6 K). Between 
p ~ 0.5 g cm~ 3 and 20 g cm -3 , Li-like ions gradually disap- 
pear with the n — 2 shell to the benefit of He-like ions. On 
the contrary (see Fig.Ub)), FLYCHK still allows recombina- 
tion into the Be-like ion (v = 4) until p ~ 8 g cm~ 3 where it 
is sharply pressure ionized, just before Li-like ions. Oxygen 
remains at higher densities in the He-like stage up to p ~ 10 3 
g cm -3 where the last Is shells for He-like and H-like ions are 
moved into the continuum. In this density range our results 
show a more complex behaviour due to several effects: the 
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FIG. 5: (Color online) Number fractions of different carbon (dotted 
lines) and oxygen (solid lines) ions in the mixture of equal number 
of ions of each chemical element at p = 1 g cm -3 as functions of T. 



FIG. 4: (Color online) (a) Number fractions of ions for a pure oxygen 
plasma (solid lines) at T — 6 x 10 J K as a function of the density 
p. Numbers of bound electrons v for each ion are marked near the 
curves. Dashed lines are twice the ionization fractions Xj V for oxy- 
gen in a plasma mixture composed of equal numbers of C and O 
nuclei at the same temperature, (b) Mean effective charge (Z) for 
the pure oxygen plasma (solid line). The dot-dashed line is the result 
from FLYCHK under the same conditions. 



free electrons become degenerate, the ionization occurs on a 
large density range due to band broadening, and the influence 
of neutrality volumes is strong. 

Dashed lines in Fig. Ua) show the ionization fractions of 
oxygen in a mixture of equal numbers of C and O nuclei. 
The fractions Xj v have been multiplied by two to be com- 
pared with the pure oxygen plasma case. The behavior of the 
ionization equilibrium of oxygen is quite similar in the two 
cases, the biggest difference occurring in the pressure ioniza- 
tion zone around p ~ 10 3 gem -3 . In Fig. [5] we show the 
temperature dependence of the partial ion fractions for the two 
chemical elements C (dotted lines) and O (solid lines) in the 
same mixture at p = 1 g cm~ 3 . At this density, the curves for 
the fraction of z^-ions for carbon are roughly shifted to lower 
temperatures when compared to the oxygen ones by the ratio 
of the ionization energies (Zc — v) 2 j (Zo — v) 2 . 

Comparing Fig. [5] and analogous figures for mixtures and 
pure substances, we notice that the temperature and density 
dependences of partial ion fractions in mixtures of different 
elements are similar to those in pure substances, scaled by the 
element abundances Yj. Two examples are presented in Fig. [6] 
Here, the dashed lines show partial fractions Xj V of oxygen 
ions with v bound electrons in a mixture where half the ions 
are oxygen ions and half are carbon ions, at (a) p = 0.1 
g cm -3 and (b) p = 10 g cm -3 , as function of T, while the 
dotted lines show fractions of the same ions, multiplied by 
0.5, in a pure oxygen plasma. The close agreement between 



these curves, which has been tested at several densities (see, 
e.g., Fig. Eta)) or compositions, prompted us to suggest the 
following method for an approximate solution of the ioniza- 
tion equilibrium problem and a construction of the EOS for 
plasma mixtures of different elements (at least in the temper- 
ature range considered in the present study): 

• calculate the fractional numbers of different ions Xj V 
for pure substances at the same n e and T (let us denote 
their values Xj V ); 

• multiply them by Yj and keep them fixed: Xj V = x® v Yf, 

• adjust the electron chemical potential p e so as to fulfill 
the last condition in Eq. (|8). 

The advantage of this method for a mixture of J chemical 
elements is that, instead of minimizing F to t in a space of 
1 + Ylj=i (Zj ~ 1) independent parameters, one needs to per- 
form J minimizations in a space of Zj parameters each, and 
then to adjust p e . As a rule, these J partial minimizations 
go much faster than the single minimization in the space of 
all parameters, thus saving the computer cost of the proce- 
dure. For example, it is much easier to find the minima of 
Ft t for carbon (6 parameters) and for oxygen (8 parameters) 
separately, than to find a minimum for their mixture (13 pa- 
rameters). Moreover, having found the number fractions of 
ions for pure substances once, one can calculate thermody- 
namic functions of mixtures with any fractions of chemical 
elements Yj without repeating the whole minimization pro- 
cedure. Note however that carbon and oxygen are chemical 
species that are not very different, so that one should proceed 
with caution before applying this simplified approach to other 
mixtures without further studies. 
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FIG. 6: (Color online) Number fractions of different oxygen ions as 
functions of T in the mixture of equal number of ions of each chem- 
ical element C and O at (a) p — 0.1 g cm -3 and (b) p = 10 g cm -3 . 
Dashed lines show the result of direct free energy minimization for 
the mixture, and dotted lines show the result of an approximation 
based on the solution of the ionization equilibrium problem for pure 
C and O plasmas (see text). 



FIG. 7: (Color online) (a) Normalized pressure PV/NiksT for pure 
carbon along the isotherm T = 5 x 10 5 K. The present results (long- 
dashed line) are compared with the results from Ref. [4] (dot-dashed 
line) and Ref. 1 1] (short-dashed line). Note the different scale in 
the figure above and below the middle horizontal line, (b) Number 
fractions of different carbon ions (present results) for the same tem- 
perature. 



B. Thermodynamic functions 

To illustrate our results for the EOS and compare them with 
the results of Ref. yl, we show in Fig.[7ja) an isotherm for the 
pressure P(p) for carbon at T — 5 x 10 5 K. Our present re- 
sult is shown by the long-dashed line, our previous result [1] 
by short dashes, and the dot-dashed line shows the result ob- 
tained with the FGV EOS [4] based on a Thomas-Fermi model 
at P ^ 3 g cm -3 . In order to make the differences visible in 
the figure, we have normalized the pressure by its value for an 
ideal Boltzmann gas of ions, Nik^T/V . The vertical scale is 
smaller for the upper part of the figure, to take into account the 
rapidly increasing pressure contribution of degenerate elec- 
trons. The differences between our and FGV results have be- 
come smaller for p > 10 g cm~ 3 than in Ref. IjJ], and the 
additional features have become less spectacular. This change 
is mainly caused by the improved treatment of the electron 
continuum in the present work (see Sec. HIT). 

However, the additional features still persist near the den- 
sities where electron shells become pressure ionized, as ex- 



pected. For instance it can be seen from the ionization frac- 
tions as given in Fig. [7jb) for the same isotherm that the fea- 
ture in the pressure just above p = 1 g cm -3 is linked to the 
disappearance of the Li-like carbon ion. These features are 
clearly revealed in the logarithmic derivative of the pressure, 
X P = (d In P/d In p)t, shown in Fig.[8ja). In this panel, the 
dot-dashed, short-dashed, and long-dashed lines correspond 
to the same three models for carbon plasma as in Fig. |7{a). 
In addition, the solid line shows Xp f° r me oxygen plasma 
at the same temperature. When analyzed with the help of 
the ionization fractions as plotted in Fig. (8jb) the features 
for carbon around p ~ 5 g cm -3 and 3 x 10 2 g cm~ 3 ap- 
pear to be linked to the pressure ionization in H-like ions of 
the n = 2 and n = 1 shells respectively; those for oxygen 
around p ~ 10 g cm~ 3 and 6 x 10 2 g cm~ 3 are linked to 
the same shells with a supplementary contribution from the 
He-like ions. In the latter case, the features due to pressure 
ionization are more pronounced because of the higher binding 
energies compared to carbon, and they are shifted to higher 
densities because of the higher net charge and respectively 
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FIG. 8: (Color online) (a) Logarithmic pressure derivative \p = 
dlnP/dlnp along the isotherm T — 3.16 x 10 6 K. The present 
data for carbon (long-dashed line) are compared with the results from 
Ref. 0] (dot-dashed line) and Ref. [ 1] (short-dashed line). The solid 
line shows the isotherm of \p f° r oxygen, (b) Present results for 
the number fractions of carbon ions (dashed lines) and oxygen ions 
(solid lines) at the same temperature. 



smaller neutrality volumes of ions (this shift can be roughly 
estimated as the ratio of the volumes of hydrogen-like carbon 
and oxygen ions, ~ (Z /Z c ) 3 = 64/27). Note also that the 
pressure ionization of a given (sub-)shell leaves an imprint on 
thermodynamics functions at different densities depending on 
which ionization stage is dominant. For instance the features 
due to the n = 2 shell in carbon appear at p ~ 1 g cm -3 
or 5 g cm -3 depending on whether Li-like or H-like ions are 
concerned (compare Fig. [7] and Fig. |H). 

Figure [9] shows three isotherms for the internal energy per 
unit mass for carbon (long-dashed lines), oxygen (solid lines), 
and the mixture containing equal numbers of C and O nu- 
clei (short-dashed and dotted lines). The energy per unit 
mass U* = U/(NiJ2j m jYj) + Uq is measured from the 
ground state energy of nonionized atoms, which corresponds 
to a shift with respect to the electron continuum level equal to 

U£ = J2 3 Y J U o,r where U lj = 8 - 28 x 1013 er § for car " 
bon and Uq j = 1.23 x 10 14 erg g _1 for oxygen. For the mix- 
ture, the short-dashed lines portray the results of the numeri- 
cal minimization of the complete free energy i*tot, while the 



14.5 = 



00 
00 



00 

o 



14 



13.5 



I I I I III 1 1 I I I I I II 1 1 I I I I III 1 f~ 



6.1 



/J, 



5.7 

C 


C + O (num.) 
C + (appr.) 

_miJ 



0.01 



0.1 



1 10 

p (g cm" 3 ) 



100 



1000 



FIG. 9: (Color online) Specific internal energy logarithms (in cgs 
units) for various isotherms: T = 5 x 10 5 K, 1.26 x 10 6 K, and 
3.16 x 10 6 K(the curves are marked by log 10 T values) for oxygen 
(solid lines), carbon (long-dashed lines), and the mixture of equal 
number of ions of each chemical element (short-dashed and dotted 
lines). In the latter case short-dashed lines show the result of direct 
free energy minimization for the mixture, while dotted lines show 
the result of an approximation based on the solution of the ioniza- 
tion equilibrium problem for pure C and O plasmas (the difference 
between them is almost unnoticeable; see text). 



dotted lines illustrate the results of the approximate method 
described in Sec. IV Al The near coincidence of the dotted and 
short-dashed lines, indistinguishable at the scale of the figure, 
confirms the accuracy of the approximate method in that case. 

As an example of a second-order thermodynamic function 
calculated with our model, Fig. [I0]presents isochores for the 
heat capacity at (a) p = 0.1 g cm~ 3 , (b) p = 1 g cm~ 3 , 
and (c) p = 10 g cm -3 for carbon (long-dashed line), oxygen 
(solid line), and the carbon-oxygen mixture with Yj = 0.5. 
For the mixture, the short-dashed lines portray the results of 
the numerical minimization of the complete free energy F to t, 
while the dotted lines illustrate the results of the approximate 
method described in Sec. IV Al As well as in Fig.|£l the short- 
dashed and dotted lines nearly coincide. For comparison, lines 
drawn with alternating short and long dashes show the results 
obtained for the fully ionized plasma model. The presence 
of bound states yields big bumps on the curves. Comparison 
with Figs.[3]and|6](for 0.1 — 0.2 g cm -3 and 10 g cm -3 ) and 
with Fig. [5] (for 1 g cm -3 ) reveals that these bumps are as- 
sociated with the thermal ionization of heliumlike ions into 
hydrogenlike ones and then into nuclei. At higher tempera- 
tures, the results for fully ionized plasmas are recovered, as it 
should be due to the complete ionization of C and O at such 
high temperatures, in agreement with the Saha equation. 
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FIG. 10: (Color online) Normalized specific heat per ion for the 
isochores (a) p = 0.1 g cm -3 , (b) p = 1 g cm" 3 , and (c) 
p = 10 g cm" 3 for carbon (long-dashed line), oxygen (solid line), 
and mixture of equal number of ions of each chemical element (short- 
dashed and dotted lines). In the latter case, the short-dashed line 
shows the result of direct free energy minimization for the mixture, 
while the dotted line shows the result of an approximation based on 
the solution of the ionization equilibrium problem for pure C and O 
plasmas. For comparison, alternating short and long dashes show 
analogous isochores for the fully ionized plasma model. 



VI. CONCLUSIONS 

We have further improved our model JTt] for the calculation 
of the EOS for dense, partially ionized plasmas, based on the 
free energy minimization method, suitable to handle pressure 
ionization regimes. The free energy model is constructed in 
the framework of the chemical picture of plasmas and includes 
a detailed, self-consistent treatment of the quantum states of 
partially ionized atoms in the plasma environment as well as a 
quantum treatment of the continuum electrons. The improve- 
ments with respect to our former pure carbon calculations in- 
clude implementing updated fitting formulae ll25l l26Tl for the 
long-range Coulomb contribution to the free energy, the use 
of the OP database for energies of bound state configurations 
dlH, an increased range of density and the extension to mix- 
tures of different chemical elements. 

We have extended the calculations from pure carbon to ar- 
bitrary carbon-oxygen mixtures. We have also suggested an 
efficient, approximate but accurate method of EOS calcula- 
tion in the case of mixtures of different chemical elements, 
based on the detailed information about the various ionic state 
fractions for the pure species. 

The present EOS results can be used in studies of carbon- 
oxygen plasmas in the domain of so-called warm dense mat- 
ter as well as in astrophysical calculations of stellar structure 
and evolution. The possible astrophysical applications include 
oxygen and carbon plasmas and carbon-oxygen mixtures in 
various types of stars, for example the interiors of carbon- 
oxygen white dwarfs. However, in order to apply it also to 
the outer parts of the white dwarf envelopes, it is desirable to 
extend the results to lower temperatures. We are planning to 
perform such extensions in the future work. 
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